package com.numericalmethod.algoquant.model.ralph2009.draft;
import java.io.BufferedReader;
import java.io.FileReader;
import java.util.Vector;

import rcaller.RCaller;
import rcaller.RCode;


public class Estimate_param {

	
	 
	String R_path = "C:\\Program Files\\R\\R-2.15.2\\bin\\x64\\Rscript.exe";
	String InputFileCSV = "stock_n_yield.csv";
	
	Calibrated_param[] cali_parm; 
	Vector<Double> vec1;
	Vector<Double> vec2;
    double[][] cov;
    double[][] timeLag1_coeff;
    double[]const_coeff;
	
	
	int iWindowSize;
	double x1[];
	double x2[];
	int iOffset;
	double del_t;
		
	void Init_All() {
		
		System.out.println("================================");
		System.out.println("Warning! The default R installatoni path is this:"+R_path);
		System.out.println("If this doesn't match your system, please modify it inside Estimate_param.java class");
		System.out.println("Also, user input is assumed to be "+InputFileCSV);
		System.out.println("Please modify accordingly.");
		System.out.println("================================");
			
		vec1 = new Vector<Double>();
		vec2 = new Vector<Double>();
						  
		x1 = new double[iWindowSize];
		x2 = new double[iWindowSize];
		
        cov = new double[2][2];
        timeLag1_coeff = new double[2][2];
        const_coeff = new double[2];		
        
        //del_t = 1.0/12.0;
        del_t = 1.0; //????  		
	}
	
	// with user input rolling window size 
	Estimate_param(int iInput_WindowSize) {
		
		iWindowSize = iInput_WindowSize;
		Init_All();
	}
	
	Estimate_param() {
		
		iWindowSize = 109;// please change this number to something reasonable.  This is just for testing.
		
		Init_All();
	}
	
	void DoEstimation() {
		ReadInputFile();
			
		
		cali_parm = new Calibrated_param[vec1.size()];
		
		for (int ii=0; ii < vec1.size(); ii++) {
			cali_parm[ii] = new Calibrated_param();
		}

		
		if ((vec1.size() - 1) >= (2*iWindowSize-1)) {
			//iOffset = (vec1.size() - 1) - (2*(iWindowSize)-1);
			iOffset = (vec1.size() - 1) - iWindowSize - (iWindowSize-1); // 1st starting pt for storing the fitting
		}
		else {
			System.out.println("Not enough data for window size= "+iWindowSize + "for data point = " + vec1.size());
			iWindowSize = vec1.size()/2; // (2*win_size-1) <= total_element-1
			System.out.println("Resized window  =" +iWindowSize);
			iOffset = (vec1.size() - 1) - iWindowSize - (iWindowSize-1); // 1st starting pt for storing the fitting
		}
				
		// do it from T-n time back to time 0.
		for (int ii=iOffset; ii >= 0; ii--) {
			CreateOneRollingWindow(ii);
			CallRToSolve(x1,x2);
		    ParseCov(cov);	        
		    ParseCoeff(timeLag1_coeff, const_coeff);
		    EsimateNow(timeLag1_coeff, const_coeff, cov, del_t,ii);
		    Print(ii);		    
		    //System.out.println("iOffset="+iOffset+" ii="+ii);
		}
	}
	
	double CalM(Vector<Double> TmpVec1,
			  int iTmpOffset) {
		double rate_ret, M=0, discount;
		
		for (int ii=0; ii < iWindowSize; ii++) {
			rate_ret = (TmpVec1.elementAt(ii+iTmpOffset) - TmpVec1.elementAt(ii+iTmpOffset+1))/TmpVec1.elementAt(ii+iTmpOffset+1);
			if (ii > 0) {
				discount = Math.exp(-ii);
			}
			else {
				discount=1.0;
			}
			
			M += discount*rate_ret;
		}
		
		return M;
	}
	
	void CreateOneRollingWindow(int iTmpOffset) {
		double u_avg=0.0;
		
		// we assume no of elem in vector is at lest 1 elem greater than windowSize
		// i.e. vec.size=5, windowsSize=4; because for a window of generating 4 return, we need 5 points.
		
		for (int ii=0; ii < iWindowSize; ii++) {
			u_avg += vec2.elementAt(iTmpOffset+ii); // sum up all first
		}
		u_avg = u_avg / iWindowSize;
		
		for (int ii=0; ii < iWindowSize; ii++) {
			//x1[ii] = vec1.elementAt(iOffset+ii);
			x1[ii] = CalM(vec1,iTmpOffset+ii);
			x2[ii] = vec2.elementAt(iTmpOffset+ii) - u_avg; // this is X, yeild - mean(yield)
		}				
	}
	
	void CallRToSolve(double []x1, 
			double []x2) {
		try {
			  // part to call R caller
			RCaller caller = new RCaller();
			RCode code = new RCode();
		    code.clear();	      	       
		      
	        
		    caller.setRscriptExecutable(R_path);		      
	        caller.cleanRCode();
	                    
	        code.clear();          
	        code.addDoubleArray("data1", x1);
	        code.addDoubleArray("data2", x2);
	                      
	        code.addRCode("library(vars)");
	        code.addRCode("y=cbind(data1,data2)");
	        code.addRCode("mod <- VAR(y,p=1, type = \"const\")");
	                 
	        code.addRCode("sink(\"temp_parse1.txt\")");
	        code.addRCode("summary(mod)");
	        code.addRCode("sink()");

	        code.addRCode("sink(\"temp_parse2.txt\")");
	        code.addRCode("coef(mod)");
	        code.addRCode("sink()");
	                    
	                    
	        caller.setRCode(code);
	        caller.runOnly();
			
		} catch (Exception e) {
		    System.out.println(e.toString());
		}				
	}
	
	void ReadInputFile() {
		try {
			BufferedReader CSVFile = 
			        new BufferedReader(new FileReader(InputFileCSV));
		
			String dataRow = CSVFile.readLine(); // Skip first line.
			dataRow = CSVFile.readLine(); // Read first line.
			  
			double dTmpX;
			
		
			while (dataRow != null){
			    String[] dataArray = dataRow.split(",");
			   
			   for (int ii=1; ii < dataArray.length; ii++) { 
				   
				   dTmpX = Double.parseDouble(dataArray[ii]);
				   
				   if (ii == 1) { // col 1
					   vec1.addElement(dTmpX);
				   }
				   else 
					   if (ii == 2){ // col 2 
						   vec2.addElement(dTmpX);
					   }
			  //    System.out.print(dTmpX + "\t"); 
			   }
			//   System.out.println(); // Print the data line.
			   dataRow = CSVFile.readLine(); // Read next line of data.
			}		 
			CSVFile.close();	
		} catch (Exception e) {
		    System.out.println(e.toString());
		}		
	}
	
	void EsimateNow(double [][]timeLag1_coeff, 
			double []const_coeff,
			double [][]cov,
			double delta_t,
			int iTimeOffset) {
		
		// first we compute the as cholesky decomposed
		// i.e. A = L*L' where L' is the conjugate tranpose
		//[s	0] X [s	 x1]   [s^2	  sx1]
		//[x1	x2]	 [ 0 x2] = [sx1	  x1^2+x2^2]
		// we can write the closed form for all the lower triangular parameters.
		
		cali_parm[iTimeOffset].sigma_s = Math.sqrt(cov[0][0]); // s^(0.5)
		cali_parm[iTimeOffset].sigma_x1 = cov[1][0]/cali_parm[iTimeOffset].sigma_s; // cov(1,0)/s
		cali_parm[iTimeOffset].sigma_x2 = Math.sqrt(cov[1][1] - cali_parm[iTimeOffset].sigma_x1*cali_parm[iTimeOffset].sigma_x1); // sqrt(cov(1,1)-x^2)
		
		// fitted VAR(1) becomes
		// Yt = {const} + {coeff}Y(t-1)
		// expand as
		// [Mt] = [(1-phi)u0*del_t]  + [-(1-phi)del_t+1 u1(1-phi)del_t] [M(t-1)]    Sigma*eps
		// [Xt]   [ 0             ]    [0               -alpha*del_t+1] [X(t-1)] + 
		//
		//   [a] + [c d]Y(t-1)
		//   [b]   [e f]
		double a,c,d,f;
		f = timeLag1_coeff[1][1];
		c = timeLag1_coeff[0][0];
		a = const_coeff[0];
		d = timeLag1_coeff[0][1];
		
		cali_parm[iTimeOffset].alpha = -(f-1.0)/delta_t;		
		cali_parm[iTimeOffset].phi = 1.0 - (1.0 - c)/delta_t;
		cali_parm[iTimeOffset].u0 = a/((1-cali_parm[iTimeOffset].phi)*delta_t);
		cali_parm[iTimeOffset].u1 = d/((1-cali_parm[iTimeOffset].phi)*delta_t);
		
		// no momentum calibration (i.e. mean reversion)
		cali_parm[iTimeOffset].u0_prime = a/delta_t;
		cali_parm[iTimeOffset].u1_prime = d/delta_t;
		
		cali_parm[iTimeOffset].alpha_prime = cali_parm[iTimeOffset].alpha; // same as mometum version				

		
	}
	
	void Print(int iTimeOffset) {
		System.out.println("time = "+iTimeOffset);
		System.out.println("sigma_s = " + cali_parm[iTimeOffset].sigma_s);
		System.out.println("sigma_x1 = " + cali_parm[iTimeOffset].sigma_x1);
		System.out.println("sigma_x2 = " + cali_parm[iTimeOffset].sigma_x2);
		System.out.println("u0 = " + cali_parm[iTimeOffset].u0);
		System.out.println("u1 = " + cali_parm[iTimeOffset].u1);
		System.out.println("alpha = " + cali_parm[iTimeOffset].alpha);
		System.out.println("phi = " + cali_parm[iTimeOffset].phi);
		System.out.println();
	}
	
	 void ParseCov(double [][]cov) {
		
		try {
	        // parse result
			BufferedReader ParseFile1 = 
			        new BufferedReader(new FileReader("temp_parse1.txt"));
		
			String dataRow2 = ParseFile1.readLine(); // Skip first line.
			dataRow2 = ParseFile1.readLine(); // Read first line.
			  		
			int iRow, iCol;
		
			while (dataRow2 != null){
			    String[] dataArray = dataRow2.split("\\s+");
			    
			    if (dataArray[0].equals("Covariance") == true) {		    
			   
			    	dataRow2 = ParseFile1.readLine(); // skip next line
			    	dataRow2 = ParseFile1.readLine(); // skip next line		    			    			    	
			    	dataArray = dataRow2.split("\\s+");		    	
			    	
			    	iRow = iCol = 0;
			    	// 1st row of covariance
			    	for (int ii=1; ii < dataArray.length; ii++) {
			    		cov[iRow][iCol] = Double.parseDouble(dataArray[ii]);
			    		iCol++;
			    	}
			    	
			    	dataRow2 = ParseFile1.readLine(); // skip next line		    			    			    	
			    	dataArray = dataRow2.split("\\s+");
			    	
			    	iRow = 1;
			    	iCol = 0;
			    	// 2nd row of covariance
			    	for (int ii=1; ii < dataArray.length; ii++) {			    		
			    		cov[iRow][iCol] = Double.parseDouble(dataArray[ii]);
			    		iCol++;
			    	}
			    }
			   
			   dataRow2 = ParseFile1.readLine(); // Read next line of data.
			}		 
			ParseFile1.close();
		
		
		} catch (Exception e) {
		    System.out.println(e.toString());
		}		
	}
	
	
	void ParseCoeff(double [][]timeLag1_coeff, double []const_coeff) {
		
		try {
	        // parse result
			BufferedReader ParseFile1 = 
			        new BufferedReader(new FileReader("temp_parse2.txt"));
		
			String dataRow2 = ParseFile1.readLine(); // Skip first line.
			  		
			int iRow, iCol;
		
			while (dataRow2 != null){
			    String[] dataArray = dataRow2.split("\\s+");
			    
			    if (dataArray[0].equals("$data1") == true) {		    
			   
			    	dataRow2 = ParseFile1.readLine(); // skip next line
			    	dataRow2 = ParseFile1.readLine(); // skip next line		    			    			    	
		    	
			    	dataArray = dataRow2.split("\\s+");		    	
			    	
			    	iRow = iCol = 0;
			    	timeLag1_coeff[iRow][iCol] = Double.parseDouble(dataArray[1]); // (0,0) time-lag-1 coeff
			    	iCol++;
			    	
			    	dataRow2 = ParseFile1.readLine(); // skip next line		    			    			    	
			    	dataArray = dataRow2.split("\\s+");
			    	
			    	timeLag1_coeff[iRow][iCol] = Double.parseDouble(dataArray[1]);// (0,1) time-lag-1 coeff
			    	
			    	dataRow2 = ParseFile1.readLine(); // skip next line		    			    			    	
			    	dataArray = dataRow2.split("\\s+");
			    	
			    	const_coeff[iRow] = Double.parseDouble(dataArray[1]);// row 0 const			   
			    }
			    else 
			    	if (dataArray[0].equals("$data2") == true) {
				    	dataRow2 = ParseFile1.readLine(); // skip next line
				    	dataRow2 = ParseFile1.readLine(); // skip next line		    			    			    	
			    	
				    	dataArray = dataRow2.split("\\s+");		    	
				    	
				    	iRow = 1;
				    	iCol = 0;
				    	
				    	timeLag1_coeff[iRow][iCol] = Double.parseDouble(dataArray[1]); // (0,0) time-lag-1 coeff
				    	iCol++;
				    	
				    	dataRow2 = ParseFile1.readLine(); // skip next line		    			    			    	
				    	dataArray = dataRow2.split("\\s+");
				    	
				    	timeLag1_coeff[iRow][iCol] = Double.parseDouble(dataArray[1]);// (0,1) time-lag-1 coeff
				    	
				    	dataRow2 = ParseFile1.readLine(); // skip next line		    			    			    	
				    	dataArray = dataRow2.split("\\s+");
				    	
				    	const_coeff[iRow] = Double.parseDouble(dataArray[1]);// row 0 const			   
			    	
			    	}
			   
			   dataRow2 = ParseFile1.readLine(); // Read next line of data.
			}		 
			ParseFile1.close();		
		
		} catch (Exception e) {
		    System.out.println(e.toString());
		}
	}

	
}
